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A  NON-RECURSIVE  INCOMPLETE  CHOLESKY  DECOMPOSITION  METHOD  FOR 
THE  SOLUTION  OF  LINEAR  EQUATIONS  WITH  A  SPARSE  MATRIX 

I.  Introduction 

The  algorithm  described  in  this  paper  is  applicable  to  diagonally 
dominant  matrices.  It  makes  use  of  the  fact  that  an  incomplete  Cholesky 
or  LU  decomposition  which  imposes  a  certain  sparsity  is  actually  an 
expansion  in  powers  of  off  diagonal  elements.  It  is  therefore  pointless 
to  compute  the  elements  of  the  tridiagonal  matrices  and  the  solution  to 
infinite  order  by  recursive  procedures.  It  is  self-consistent  and  suf¬ 
ficient  to  compute  these  quantities  to  the  same  order  as  the  one  intro¬ 
duced  by  the  imposed  sparsity  of  the  tridiagonal  decomposition. 

The  first  section  describes  the  method  for  symmetric  matrices.  The 
extension  to  non-symmetric  matrices  is  achieved  by  iterating  the  solu¬ 
tions  for  the  symmetric  component  of  the  matrix  where  the  source  terms 
contain  the  non-symmetric  contributions  from  the  previous  iteration. 

The  second  section  gives  the  times  and  number  of  iterations  used 
for  a  test  case  taken  from  a  NRL  electrostatic  code.  It  involves  the 
solution  of  an  elliptic  partial  differential  equation  with  variable  co¬ 
efficients  in  the  two  dimensions. 

The  third  section  gives  some  conclusions. 

II.  Algorithm 

1)  General  considerations 

The  linear  system  of  equations  to  be  solved  for  the  vector  x 
can  be  written  in  the  form 

Mx  ■  y  (1) 

x,y  are  vectors  of  length  N  and  Det(M)  f  0. 
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In  order  to  explain  the  following  concepts  one  must  introduce  the 


definition  of  weakly  and  strongly  diagonally  dominant  matrices.  A 
matrix  A  is  defined  to  be  weakly  diagonally  dominant  if 


’U 


Kjl 

“  1  “  min(|  A^j  |J  Ajj|  )  *  0 


(2) 


for  all  i,j;i  /  j 

A  matrix  A  is  defined  to  be  strongly  diagonally  dominant  if 


S j  (A)  =  1  - 


£  (ble  V***1)  /|Ali>  *  °  (3) 


for  all  j . 

Let  L  and  U  be  two  lower  and  upper  triangular  matrices  respectively, 
subject  to  the  following  two  conditions 


a)  W(L~1MU*”1)  *  W(M)  *  0 


(4) 


b)  The  inversions 


L-1y 


x  =  U 


(5) 


with  z  a  vector  of  length  N  can  be  performed  exactly. 


Equation  (1)  can  be  transformed  into 


L_1Mlf 1 (UX)  -  L-1y 


(6) 


It  is  obvious  that  Eq.  (6)  is  easier  to  solve  than  Eq.  (1).  Any  itera¬ 
tion  scheme  will  use  fewer  iterations.  One  should  also  remark  that 
neither  L,  U,  or  L  U  ^  have  to  be  known,  but  only  the  results  of 
L  U  *  applied  to  a  vector. 

The  real  question  is  how  one  finds  good  matrices  L,  U  with  a 
minimum  of  operations,  such  that  Eq.  (6)  can  be  solved  with  a  few 


Iterations.  This  is  especially  important  in  time-dependent  problems 
where  a  good  approximation  to  the  solution  is  known  from  the  previous 
timestep.  It  is  clear  that  with  more  operations  one  could  find  a  bet¬ 
ter  L  and  U,  thereby  reducing  the  number  of  iterations  needed  to  solve 
Eq.  (6).  The  problem  is  to  find  an  algorithm  which  minimizes  the  total 
number  of  operations. 

For  most  physical  problems  (elliptic  equations  in  two  or  three 
dimensions)  the  matrix  M  is  sparse  with  an  (average)  bandwidth  B  «  N. 
The  number  of  operations  Op  should  be 

Op  »  B*N  f(N)  (7) 


where  f(N)  is  a  weakly  increasing  function  of  N.  Also,  the  algorithm 
should  not  contain  recursive  formula,  which  Invoke  scalar  operations  on 
a  vector  computer  (ASC,  CRAY,  ,  ,  ,)  cr  parallel  computer  (ILLIAC). 

2)  Expansion  in  connections 

In  order  that  the  following  approximation  for  L  and  U  be  valid 
the  matrix  M  must  have  two  properties 

a)  the  average  bandwidth  B  «  N 

b)  M  must  be  expandable. 

Condition  b)  will  be  explained  now. 

A  connection  of  order  is  to  be  defined  as 
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€ij  VMJi  *  V  *  1 


(8) 


with 


eik  Mik 


(9) 
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Then  M  is  expandable  if  for  all  v 

l®^1!  *  oax(|«ijl  »lMjkl  )  (10> 

It  is  obvious  that  tne  condition  of  strong  diagonal  dominance 


S(M)  >  0 

is  sufficient.  Furthermore,  if  M  is  expandable  then 
W(M)  >  0 

M  is  weakly  diagonally  dominant.  Therefore,  the  necessary  and 
sufficient  conditions  for  expandability  lie  somewhere  between  weak  and 
strong  diagonal  dominance.  The  computation  of  L,  U  will  be  an  expansion 
in  the  number  of  connections,  thus  avoiding  recursion  procedures.  In 
general,  given  a  recursion  formula 


°i  + 


ep 


i-1 


i  -  l,n 


with  o  <  €  <  1. 

The  iterative  formula 

-  dj  +  en 


v-1 

i-1 


i  *  l,n 


o 


(ID 


(12) 


will  give  the  same  results  after  n  iterations.  Terminating  after  the 
iteration  gives  the  results  to  order  eV+^.  All  recursion  formula 
which  do  occur  will  be  replaced  by  an  iteration  with  a  very  small  v. 

The  matrices  L,  U  will  be  computed  to  a  certain  order  in  the  expansion  in 
the  number  of  connections.  In  order  to  explain  this  procedure  the  tri¬ 
angular  decomposition  is  briefly  described  here  without  proof. 
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.rr  *i'  p>  *S  **MP 


Let 


1/Dj  "  Lii  “  Uil  "  Mii  ~  LikDfcUfci 


Lii  *  0  i  <  i 


Lii  *  Mji  "  ^  Ljk°kUki  j  *  1 

V  “  0  j  >  1 

v  -  Mji  -  L  Vkuki  j  *  1* 


Then 

"u  -pA1.) 

and  one  can  solve 

Z1  *  Di*yi  “  LikZk^ 

Xi  “  *i  '  UikVk 


(14) 


(15) 


directly. 

The  incomplete  Cholesky  or  incomplete  LU  decomposition  consists  of 
imposing  a  specified  sparsity  on  L,  U.  The  order  of  approximation  is 


n  =  v(s)  +  1  (16) 

where  v^  Is  the  highest  possible  connectivity  given  by  the  sparsity 
Imposed  on  L  and  U  (or  n  is  first  one  neglected). 
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The  recursion  formulae  in  computing  L,  U  (Eq.  (13))  are  reset  by 


iteration  to  the  order  ,  which  leaves  an  error  of  0(p.  +  1). 
iterations  can  be  written  as  follows 


These 


■  v  -  ,s  c  -r1  <il 


k<i 


with 


* 


V-1  DV_1  UV_1 
ji  ^  “Jk  °k  Uki 


L-  “  M,.  -  X)  L 


v  =  l.u. 


(17) 


with 


TTV  _  M  v  Tv_1  nv_1  TTV_1 

V  '  V  ‘  J4  ljk  \  "kl 


1/D°  =  M. 


ii 


*  Mji 

j 

is  i 

“ii 

=  Mji 

j 

£  i 

V 

zi 

*  D?  (yA  - 

L 

k<i 

L?f  z; 
ik  1 

V 

xi 

-■S-E 

k>i 

Uik 

U  V' 

Dk  *k 

V  =  l,n 


(18) 


6 i 
o 


■S  yi 


(19) 


This  procedure  gives  a  matrix 
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whose  off  diagonal  elements  consists  of  errors  of  order  +  1 

left  by  imposing  a  specified  sparsity  S,  and  errors  of  order  t*  +  1 
caused  by  resetting  the  recursion  formula  by  iteration. 

3)  Example 

As  an  example  the  five  point  difference  formula  for  an 
elliptic  equation  in  two  dimensions  is  taken.  M  at  one  point,  j,k  is 
then 


M  =  -  u 

X 

-  U. 

l 

"  ‘Sc 

k-1 

1  k 

k+1 

j-l 

j 

j+l 

imposing  the 

same  sparsity  on  ' 

L  and  U  a 

L(V>  -  0 

1/D<V) 

-  ‘Sc 

k-1 

k 

-  u 

y 

k+1 

j-l 

j 

j+1 

““y 

k-1 

X 

3. 

1 

It 

> 

i/d<v> 

0 

k 

0 

k+1 

k-1 

k 

k+1 

*  l 

d<2>  - 

1 

1-p  2-p,  2 

rx  Ky 

(20) 


(21) 


(22) 


(23) 
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This  yields: 


This  is  not  a  good  approximation  to  a  unity  matrix.  Therefore,  one  more 
iteration  step  has  to  be  taken! 


J-3  J-2  J“1  5  J+l  j+2  j+3 


This  is  a  much  better  approximation  to  the  unity  matrix  correct  to 

terms  of  third  order,  except  for  the  off  diagonal  terms  p  p  caused 

x  y 

By  the  sparsity  imposed  on  L  and  U,  A  complete 


8 


The  next  iteration  vectors  are  then  given  by 


(x+1) 


x<»>  +  a^.p(» 

B(x>  * 


.(X+l)  .  „(X)  .  «£i  „(X) 
6(X)  " 


,<X+1)  .  (A+l)  .  »(H'1)  „(X) 

'  "  S  (X)  p 


(30) 


5)  Non-symmetric  matrices 

Non-symmetric  matrices  M  arise  in  physical  problems  from 
different  sources:  gradients  that  can  also  be  present  in  elliptic 
equations,  non-separable  coordinate  systems  and  most  important  from 
Neumann  boundary  conditions.  The  method  proposed  by 

Kershaw  did  not  work  very  well  in  test  problems.  It  consists  essen¬ 
tially  of  multiplying 

L  1MU~1(Ux)  =  L-1y  (31) 

-1  -1  ^ 

by  the  transposed  matrix  (L  MU  )  and  then  solving  the  resulting  linear 

system.  One  can  easily  see  that  off  diagonal  elements  are  multiplied 

essentially  by  two.  Therefore,  the  convergence  of  the  conjugate 

T 

gradient  method  is  slowed.  Even  worse  is  to  take  M  M.  Then  the  con¬ 
dition  of  expandability  is  not  fullfilled  for  elliptic  equations.  The 
method  used  here  is  to  solve  for  the  asymmetry  by  iteration.  Let 


M,  -  ±  (M  +  MT) 

6M  -  ^  (M  -  MT) 

Let  a  denote  the  iteration  for  the  asymmetry,  then 

M  x(0)  ■  y  -  6Mx^°  ^ 
s 


(32) 


(33) 


i 


4 


10 


One  can  go  one  step  further  and  correct  x^°  ^  as  used  in  the  above 
equation.  Let 


Then 


or 


t-(0>  -  S  -  Mx(°>  -  MC*0-1) 


Mg(ra  +  6x)  *  0 


M  6x  =  - 
s 


r 


(a) 


(34) 


(35) 


Now  use  LL  in  above  equation  and  define 

-  (LLT)  r 


(a)  =  „(<>)  _  /ttTn_1_(o) 


for  the  right  hand  side  of  Eq.  (33). 


(36) 


III.  Numerical  results 

Test  runs  have  been  made  with  an  elliptic  equation  which  arises 
in  numerical  simulations  of  electrostatic  plasma  flow: 

(37) 


/  a_  a  .  ,  a  \  f  _  So 

\3x  dx  &)'  J  5y  I  dx 

'  '  0  <  x,y  <  1 

[-(=?)'•«] 


cf  =  1  +  g*exp 


The  partial  differential  equation  was  translated  into  a  five  point 
formula  in  two  different  ways;  a)  by  leaving  o  inside  of  the 
second  derivative  and  b)  by  taking  the  differentiation  of  o  out  of  the 
second  derivatives  and  treating  the  first  derivatives  of  o  separately. 
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pi*ULRI*PI 


while 


dx?  [(aj+l,k  +  aj,k)  (fj+l,k  "  ‘  (aj,k 

dy2  [(aj,k+l  +  (fj,k+l  ‘  fj .k*  "  (oj,k 

2dx  (CTj+l,k  "  aj-l,k^ 
b)  gives 

d^2  [(fj+l,k“2fj.k  +  fj-l,k)+  4^  (°j+l,k  ' 
dy2  [(fj,k+l~2fj,k  +  fj,k-l)  +  4^k  (aj,k+l 

2a^kdx  (°j+l,k  “  ^-l.k^ 


+  ^-1,^  (fj,k'  f j-l.k^J 

+  °j,k-l)  (fj,k"  fj,k-l>] 
(38) 


°j-l,k)(fj+l,k  “  Vi.d 

0j,k-l)(fj,k+l"  f j ,k-l)] 
(39) 


Version  a)  is  symmetric;  b)  is  not. 

4 

Test  runs  have  been  made  with  8  from  10  to  .i,  SL  ,  i.  from  .1  to 

x  y 

.5.  The  number  of  points  in  x  and  y  has  been  varied  between  25 
and  100.  The  rras  error: 

rms  5  varied  between  10  ^  and  10  ^ . 

s2 

The  general  experience  gained  by  the  test  runs  can  be  summarized  as 
follows. 

1)  The  need  for  double  precision  on  the  ASC  (the  ASC  has  a 
32-bit  word).  The  reason  seems  to  be  that  the  orthogonali- 
zatlon  has  to  be  achieved  with  high  precision. 

2)  A  simple-minded  iterative  scheme  with  single  precision  worked 
well  only  for  a  low  accuracy  and  relatively  small  n  ^  30. 

The  scheme  eventually  did  converge  but  with  a  great  number  of 


Iterations. 
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3)  The  use  of  recursion  rather  than  iteration  for  the  incomplete  LU 
decomposition  for  a  specific  approximation  did  not  (essentially) 
change  the  convergence  rate. 

4)  The  use  of  n  =  3,  =  2  compared  to  p.  ■  2,  =  1  decreases 

the  number  of  iterations  in  about  the  same  ratio  as  the  number  of 
operations  per  point  increases.  Therefore,  the  total  time 
remained  essentially  constant. 

5)  The  use  of  recursion  in  the  second  index  (allowing  partial 
vector ization  on  a  vector-computer)  only  decreased  slightly  the 
number  of  iterations.  It  seems  to  be  more  appropriate  -  at  least 
in  the  test  cases  run  -  to  use  the  same  approximation  in  the  x  and 
y  directions.  There  may  be  asymmetric  cases  where  this  will  not 
be  true. 


Timing 


The  number  of  operations  per  point  is: 

Scalar  products 
Compute  x,r,p 
q  =  Mp 


Total 


rms  error 


The  number  of  operations  for  x  ■  (LL1)  depends  on  the  approximation. 
Let  nsca  be  the  equivalent  number  of  vector  operations  for  one  scalar 
operation;  then  the  ratio  of  operations  for  tne  equivalent  number  for  a 
recursive  procedure  becomes  for  p.  ■  2,  "  1»  for  ■  3,  Vg  ■  2 


1  4nsca  +  24 


L2  4nsca  +  26 
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T 


T 


This  gives  for  nsca  *  20  (»  factor  for  the  ASC) 
rL  -  .37  r2  -  .55 


The  iterative  procedure  can  be  vectorized  over  the  whole  array.  For  the 
recursive  procedure  the  remaining  vectorization  can  be  only  achieved  for 
an  inner  loop,  thus  increasing  the  setup  times  compared  to  the  iterative 
procedure.  Also  the  more  efficient  use  of  two  pipes  on  the  ASC  machine 
decreases  the  ratio  further.  Test  runs  have  shown  for  n  ■  2,  v  ■  1  an 
overall  savings  of  about  a  factor  of  5. 

The  total  number  of  iterations  seems  to  be  proportional  to 

N  «  (N  N  )3/2 
i  x  y 


The  dependence  of  on  the  error  reduction  rate  seems  to  be  more  compli¬ 
cated,  roughly  speaking  5  iterations  per  factor  10,  such  that  the 
convergence  factor  is  fre<j  *  (.1)^^.  In  contrast  to  ADI  and  other 
methods  the  convergence  factor  seems  to  be  more  or  less  constant  and 
independent  of  the  error  itself. 


Imposing  the  Neumann  boundary  condition  ^  *  0  at  y  *  ±  1 
introduces  an  asymmetric  matrix  M.  The  number  of  iterations  about  the 
asymmetry  is  on  average  two.  The  program  imposes  at  each  iteration  an 
error  limit  which  is  about  1/2  the  error  of  the  asymmetry.  It  thereby 
avoids  unnecessary  iterations  for  the  symmetric  solutions.  Test  runs 
have  shown  that  the  number  of  iterations  (computing  time)  increases  by 
about  50%,  when  compared  to  the  same  problem  using  Dirichlet  boundary 
conditions . 

In  time  dependent  problems  the  computing  time  depends  to  a  large 
extent  on  the  guess  of  fQ.  Crude  time  dependent  calculations  where  the 
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center  of  the  exponential  function  is  simply  shifted  and  f  is  set  to 
the  previous  solution,  gives  a  reduction  of  about  3-10  over  that  given 
by  f  ■  0,  depending  on  the  required  rms  error  and  Nx,  Ny.  The  relation¬ 
ships  for  and  f  ^  still  hold  approximately.  The  reduction  is  about 
10  for  higher  accuracy  and  lower  for  a  lower  accuracy  because  the 
change  in  f  is  the  same. 

S.  Zalesak  has  used  the  scheme  extensively  in  his  electrostatic  code. 

_3 

For  small  error  (rms  *  10  )  and  bad  approximations  the  code  runs  about 

-4 

as  fast  as  a  vectorized  ADI.  For  higher  accuracy  (10  )  and  a  good 

approximation  the  code  runs  about  two  times  faster  than  ADI  with  an 

_3 

accuracy  of  10  .  ADI  did  not  converge  after  a  reasonable  number  of 

-4 

iterations  for  a  rms  error  of  10 

IV.  Conclusion 

A  vectorized  incomplete  Cholesky  description  scheme  for  the 
solution  of  large  linear  systems  for  sparse  matrices  has  been  developed. 
The  vectorization  is  achieved  by  systematically  replacing  recursive 
formulae  by  non-recursive  expansions  in  connection  strength,  both  in  the 
incomplete  LU  decomposition  and  in  applying  (LU)  *  to  a  vector .  The  con¬ 
jugate  gradient  method  assures  convergence.  The  saving  in  computing  time 
over  the  recursive  solution  on  a  vector  machine  is  about  a  factor  five. 
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OICY  ATTN  G.  HYDE 


CORMCLL  UNIVERSITY 

fsasFE  ^485octr,cal 


OICY  ATTN  D.  T.  FARLEY  JR 
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electrospa.ce  systems,  inc. 

SOX  1359 

RICHARDSON,  Tx  75080 

01CY  ATTN  h.  LOGSTON 

01CY  ATTN  SECURITY  (PAUL  PHILLIPS) 

ESL  INC. 

495  JAVA  DRIVE 
SUNNYVALE,  CA  94086 

01CY  ATTN  J.  ROBERTS 
OICY  ATTN  JAMES  MARSHALL 
01CY  ATTN  C.  W.  PRETTIE 

PORO  AEROSPACE  i  COMMUNICATIONS  CORP 
3939  P ASIAN  WAY 
PALO  ALTO,  CA  94505 

OICY  ATTN  J.  T.  MATTINGLEY 

GENERAL  ELECTRIC  COMPANY 
SPACE  OIVISION 
VALLEY  FORGE  SPACE  CENTER 
GODDARD  BLVO  KING  OF  PRUSSIA 
P-  0.  BOX  8555 
PHILADELPHIA,  PA  19101 

OICY  ATTN  M.  H.  BORTNER  SPACE  SCI  LAS 

GENERAL  ELECTRIC  COMPANY 
P.  0.  BOX  1122 
SYRACUSE,  NY  15201 

OICY  ATTN  F.  RE  1  BERT 

GENERAL  ELECTRIC  COMPANY 
TEMPO-CENTER  FOR  ADVANCED  STUDIES 
816  STATE  STREET  (P.O.  DRAWER  QQ) 

SANTA  BARBARA,  CA  95102 
OICY  ATTN  DASIAC 
OICY  ATTN  OON  CHANDLER 
OICY  ATTN  TOM  BARRETT 
OICY  ATTN  TIM  STEPHANS 
OICY  ATTN  WARREN  S.  KNAPP 
OICY  ATTN  WILLIAM  MCNAMARA 
OICY  ATTN  B.  GAMS ILL 
OICY  ATTN  MACK  STANTON 

GENERAL  ELECTRIC  TECH  SERVICES  CO.,  INC. 

HMES 

COURT  STREET 
SYRACUSE,  NY  15201 

OICY  ATTN  G.  MILLMAN 


GENERAL  RESEARCH  CORPORATION 
SANTA  BARBARA  DIVISION 
P.  0.  BOX  6770 
SANTA  BARBARA,  CA  95111 

OICY  ATTN  JOHN  ISE  JR 
OICY  ATTN  JOEL  GARBARINO 

GEOPHYSICAL  INSTITUTE 
university  of  Alaska 
FAIRBANKS,  AK  99701 

(ALL  CLASS  ATTN:  SECURITY  OFFICER) 
OICY  ATTN  T.  N.  DAVIS  (UNCL  ONLY ) 

01CT  ATTN  NEAL  BROWN  (UNCL  ONLY) 

OICY  ATTN  TECHNICAL  LIBRARY 

GTE  STL VANIA,  INC. 

ELECTRONICS  SYSTEMS  GRP-EASTERN  DIV 
77  A  STREET 
NEEOHAM,  MA  02194 

01CT  ATTN  MARSHAL  CROSS 


ILLINOIS,  UNIVERSITY  OF 
DEPARTMENT  OF  ELECTRICAL  ENGINEERING 
URBANA,  1L  61805 

OICY  ATTN  K.  YEH 

ILLINOIS,  UNIVERSITY  OF 
107  COBLE  HALL 
SOI  S.  WRIGHT  STREET 
URBANA,  IL  60680 

CALL  CORRES  ATTN  SECURITY  SUPERVISOR  FOR) 
OICY  ATTN  K.  YEH 

INSTITUTE  FOR  DEFENSE  ANALYSES 
400  ARMY-NAVY  DRIVE 
ARLINGTON,  VA  22202 

OICY  ATTN  J,  M.  AE IN 
OICY  ATTN  ERFCST  BAUER 
OICY  ATTN  HANS  WOLFHARD 
OICY  ATTN  JOEL  BENGSTON 

HSS,  INC. 

2  ALFRED  CIRCLE 
BEDFORD,  MA  01750 

OICY  ATTN  DONALD  HANSEN 

INTL  TEL  t  TELEGRAPH  CORPORATION 
500  WASHINGTON  AVENUE 
NUTLET,  NJ  07110 

OICY  ATTN  TECHNICAL  LIBRARY 

JAYCOR 

1401  CAMINO  DEL  MAR 
DEL  MAR,  CA  92014 


JOHNS  HOPKINS  UNIVERSITY 
APPLIED  PHYSICS  LABORATORY 
JOHNS  HOPKINS  ROAD 
LAUREL,  TC  20810 

OICY  ATTN  DOCUMENT  LIBRARIAN 

oicy  attn  Thomas  potemra 

OICY  ATTN  JOHN  DASSOULAS 

LOCKHEED  MISSILES  t  SPACE  CO  INC 
P.  0.  BOX  504 
SUNNYVALE,  CA  94088 

OICY  ATTN  DEPT  60-12 
OICY  ATTN  D.  R.  CHURCHILL 

LOCKHEED  MISSILES  AND  SPACE  CO  INC 
5251  HANOVER  STREET 
PALO  ALTO,  CA  94304 

OICY  ATTN  MARTIN  WALT  DEPT  52-10 
OICY  ATTN  RICHARD  G.  JOHNSON  DEPT  52-12 
OICY  ATTN  W.  L.  IMMOF  DEPT  52-12 

KAMAN  SCIENCES  CORP 
P.  0.  BOX  7463 
COLORADO  SPRINGS,  CO  80933 
OICT  ATTN  T.  MEAGHER 

LINCABIT  CORP 
10453  ROSELLE 
SAN  DIEGO,  CA  92121 

OICT  ATTN  IRWIN  JACOBS 

M.I.T.  LINCOLN  LABORATORY 
P.  0.  BOX  73 
LEXINGTON,  MA  02173 

OICT  ATTN  DAVID  M.  TOWLE 

OICT  ATTN  L.  LOUGHLIN 


MARTIN  MARIETTA  CORP 
ORLANDO  DIVISION 
P.  0.  BOX  5837 
ORLANDO,  FL  32805 

OICY  ATTN  R.  HEFFNER 


MCOOANELL  DOUGLAS  CORPORATION 
5301  BOLSA  AVENUE 
HUNTINGTON  BEACH,  CA  92647 
01CY  ATTN  N.  HARRIS 
01CY  ATTN  J.  MOULE 
01CY  ATTN  GEORGE  MROZ 
01CY  ATTN  W.  OLSON 
01CY  ATTN  R.  W.  HALPRIN 
01CY  ATTN  TECHNICAL  LIBRARY  SERVICES 


MISSION  RESEARCH  CORPORATION 
73S  STATE  STREET 
SANTA  BARBARA,  CA  93101 
01CY  ATTN  P.  FISCHER 
01CY  ATTN  w.  F.  CREVIER 
01CY  ATTN  STEVEN  L.  GUTSCHE 
01CY  ATTN  D.  SAPPED  I  ELD 
01CY  ATTN  R.  BOGUSCH 
01CY  ATTN  R.  HEMIR1CK 
01CY  ATTN  RALPH  <IL8 
01CY  ATTN  DAVE  SOWLE 
01CY  ATTN  F.  FAJEN 
01CY  ATTN  M.  SCHEI8E 
01CY  ATTN  CONRAD  L.  LONGMIRE 
OICY  ATTN  WARREN  A.  SCHLUETER 

MITRE  CORPORATION,  THE 
P.  0.  BOX  208 
BEDFORD,  MS  01730 

OICY  ATTN  JOMV  MORGANSTERN 
OICY  ATTN  G.  HARDING 
OICY  ATTN  C.  E.  CALLAHAN 

MITRE  CORP 

WESTGATE  RESEARCH  PARK 
1820  DOLLY  MADISON  BLVD 
MCLEAN,  VA  22101 

OICY  ATTN  W.  HALL 
OICY  ATTN  W.  FOSTER 

PACIFIC-SIERRA  RESEARCH  CORP 
1456  CLOVERFIELD  BLVO. 

SANTA  MONICA,  CA  90404 

OICY  ATTN  E.  C.  FIELD  JR 

PENNSYLVANIA  STATE  UNIVERSITY 
IONOSPHERE  RESEARCH  LAB 
318  ELECTRICAL  ENGINEERING  EAST 
UNIVERSITY  PARK,  PA  16802 

(NO  CLASSIFIED  TO  THIS  ADDRESS) 

OICY  ATTN  IONOSPHERIC  RESEARCH  LAB 

PHOTOPETR1CS,  INC. 

442  MARRETT  ROAD 
LEXINGTON,  MA  02173 

OICY  ATTN  IRVING  L.  KOFSKY 

PHYSICAL  DYNAMICS  INC. 

P.  0.  BOX  5027 
BELLEVUE,  WA  98009 

OICY  ATTN  E.  J.  FREMOUW 

PHYSICAL  DYNAMICS  INC. 

P.  0.  BOX  1069 
BERKELEY,  CA  94701 

OICY  ATTN  A.  THOMSON 


RID  ASSOCIATES 

P.  0.  BOX  9695 

MARINA  DEL  REY,  CA  90291 

OICY  ATTN  FORREST  GILMORE 
OICY  ATTN  BRYAN  GA8BARD 
OICY  ATTN  WILLIAM  8.  WRIGHT  JR 
OiCY  ATTN  ROBERT  F.  LELEVIER 
OICY  ATTN  WILLIAM  J.  KARZAS 
OICY  ATTN  H.  ORY 
OICY  ATTN  C.  MACDONALD 
OICY  ATTN  R.  TURCO 


RAM)  CORPORATION,  THE 
1700  MAIN  STREET 
SANTA  MONICA,  CA  90406 

OICY  ATTN  CULLEN  CRAIN 
OICY  ATTN  ED  BEDROZ I AN 

RIVERSIDE  RESEARCH  INSTITUTE 
80  WEST  EM)  AVENUE 
NEW  YORK,  NY  10023 

OICY  ATTN  VINCE  TRAPANI 

SCIENCE  APPLICATIONS,  INC. 

P.  0.  BOX  2351 
LA  JOLLA,  CA  92038 

OICY  ATTN  LEWIS  M.  LINSON 
OICY  ATTN  DANIEL  A.  HAMLIN 
OICY  ATTN  D.  SACHS 
OICY  ATTN  E.  A.  STRAKER 
OICY  ATTN  CURTIS  A.  SMITH 
OICY  ATTN  JACK  MCDOUGALL 

RAYTHEON  CO. 

528  BOSTON  POST  ROAD 
SUDBURY,  MA  01776 

OICY  ATTN  BARBARA  ADAMS 

SCIENCE  APPLICATIONS,  INC. 
HUNTSVILLE  DIVISION 
2109  W.  CLINTON  AVENUE 

suite  700 

HUNTSVILLE,  AL  35805 

OICY  ATTN  DALE  H.  DIVIS 

SCIENCE  APPLICATIONS,  INCORPORATED 
8400  WESTPARK  drive 
MCLEAN,  VA  22101 

OICY  ATTN  J.  COCKAYNE 

SCIENCE  APPLICATIONS,  INC. 

80  MISSION  DRIVE 
PLEASANTON,  CA  94566 
OiCY  ATTN  SZ 


SRI  INTERNATIONAL 
333  RAVE NS WOOD  AVENUE 
MENLO  PARK,  CA  94025 

OICY  ATTN  DONALD  NEILSON 
OICY  ATTN  ALAN  BURNS 
OICY  ATTN  G.  SMITH 
OICY  ATTN  L.  L.  COBB 
OICY  ATTN  DAVID  A.  JOHNSON 
OICY  ATTN  WALTER  G.  CHESNUT 
OICY  ATTN  CHARLES  L.  RINO 
OICY  ATTN  WALTER  JAYE 
OICY  ATTN  M.  BARON 
OICY  ATTN  RAY  L.  LEADABRAND 
OICY  ATTN  G.  CARPENTER 
OICY  ATTN  G.  PRICE 
OICY  ATTN  J.  PETERSON 
OICY  ATTN  R.  HAKE,  JR. 

OICY  ATTN  V.  GONZALES 
OICY  ATTN  D.  MCDANIEL 
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TECHNOLOGY  INTERNATIONAL  CORP 
75  WIGGINS  AVENUE 
BEDFORD,  HA  01730 

OICY  ATTN  w.  P.  BOQUIST 

TRW  DEFENSE  t  SPACE  SYS  GROUP 

ONE  SPACE  PARK 

REDONDO  BEACH,  CA  90278 

01CY  ATTN  R.  K.  PLEBUCH 
01CY  ATTN  S.  ALTSCHULER 
01CY  ATTN  D.  DEE 

VISIDYNE,  INC. 

19  THIRD  AVENUE 

north  west  industrial  park 

BURLINGTON,  HA  01803 

oicy  attn  Charles  Humphrey 
oIcy  attn  j.  w.  carpenter 
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